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Abstract 

Four distinct theories describing the fiexural motion of thermoelastic thin plates are compared. The 
theories are due to Chadwick p^, Lagnese and Lions [S], Simmonds [3] and Norris 0]. Chadwick's theory 
requires a 3D spatial equation for the temperature but is considered the most accurate as the others are 
derivable from it by different approximations. Attention is given to the damping of fiexural waves. Analytical 
and quantitative comparisons indicate that the Lagnese and Lions model with a 2D temperature equation 
captures the essential features of the thermoelastic damping, but contains systematic inaccuracies. These 
are attributable to the approximation for the first moment of the temperature used in deriving the Lagnese 
and Lions equation. Simmonds' model with an explicit formula for temperature in terms of plate deflection 
is the simplest of all but is accurate only at low frequency, where the damping is linearly proportional to the 
frequency. It is shown that the Norris model, which is almost as simple as Simmond's, is as accurate as the 
more precise but involved theory of Chadwick. 
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Introduction 

A recent renewal of interest in thermoelasticity of vibrating structures is motivated in part 
by the observation that thermoelasticity can be the dominant and unavoidable source of 
dissipation and noise in certain circumstances. Structures both small and large in size 
can be effected. For instance, carefully controlled experiments have demonstrated that the 
dominant loss in silicon micromelectromechanical (MEMS) double paddle oscillators is due 
to thermoelasticity j^j. Conversely, large masses used in devices planned for gravitational 
wave measurements are susceptible to thermoelastic effects, which need to be controlled if 
the desired sensitivity is to be obtained [Hj. The comparative advantages of several related 
theories for the the dynamics of thermoelastic thin plates are considered here, with a view 
towards improved simulations and understanding of thermoelastic damping in devices. 

Although the subject is thin plates, it should be noted that the picture for thin beams 
is well understood. Damping of flexural motion in beams caused by irreversible thermoe- 
lastic coupling is well described by Zener's [7] theory of 1938. According to Zener, the 
instantaneous or isentropic bending stiffness of the beam relaxes to the isothermal value as a 
result of diffusion of heat produced by the alternating compression and extension on opposite 
faces of the thin structure. Zener correctly surmised that diffusion is predominantly in the 
through-thickness direction, with little or no lateral (in-plane) heat conduction, although he 
did not provide a rigorous justification for this assumption. He also derived the full form of 
the beam damping in terms of an infinite set of characteristic relaxation times defined by 
the eigenvalues of the 1-dimensional heat equation across the thickness. Similar and more 
general results for the beam were obtained by Alblas jHl El; and more recently by Lifshitz 
and Roukes [TUj . 

The situation for thin plates is more complicated by virtue of the two dimensional nature, 
which introduces the possibility of flexure in two directions. The goal of researchers has 
been to derive reduced governing equations which include thermal effects but are similar 
to the classical Kirchhoff equations. The first, and to date most successful work in this 
respect, was Chadwick pQ who provided a more rigorous basis for thermoelastic equations of 
thin beams based on the underlying 3-dimensional theory. The plate equations of Chadwick 
include 3-dimensional variability in the temperature field. This level of complexity is however 
unnecessary, as will be made clear, since the temperature variations in the plane of the plate 
are insignificant in comparison with the diffusion driven temperature gradient across the 
thickness. 

A distinct theory for thermoelasticity of thin plates was developed by Lagnese and Lions 
j2]. The Lagnese and Lions theory differs from Chadwick's in that it retains a separate 
temperature equation, although this equation does not contain variations in the through 
thickness. These are eliminated by taking the first moment of temperature across the plate 
thickness, so that both the temperature and the deflection equation are two-dimensional in 
space, and coupled. These equations bear some resemblance to those of Chadwick, and the 
connection will be discussed further below. The equations of Lagnese and Lions have been 
widely used as the starting point for studies of thermoelastic plates, e.g. [TTllT^ITBllTlllTollTol 
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Many of these works are not concerned with the mechanical applications per se but are 
focused on mathematical control issues that arise from the unique form of the coupled system 
of partial differential equations. Some, e.g. [13], have derived similar sets of equations based 
on the methods used by Lagnese and Lions j2] . 

A third theory is due to Simmonds [Sj, who derived a single equation governing the plate 
deflection. His equation contains a term that acts like a viscous damper but encapsulates 
the thermoelasticity explicitly. Simmonds proposed his equation as a simpler alternative to 
the coupled equations of Lagnese and Lions, although he did not compare his model with 
Chadwick's. 

The fourth thin plate theory considered is the most recent, due to Norris [3]. It is also 
a single equation for the plate deflection, slightly more complicated than Simmonds', as we 
will see. The Norris equation is based on recent work by Norris and Photiadis [H] who 
considered thermoelastic damping in arbitrary structures. They derived Zener-like results 
based on the asymptotically small thermal coupling. They also applied their general theory 
to time harmonic oscillation of thin plates and derived a modified plate equation without 
temperature appearing explicitly; instead it enters via a frequency dependent viscoelastic 
term. 

The purpose of this paper is to explore the similarities and differences between these 
models for thermoelastic flexural motion in thin plates. No such study has appeared to date, 
and the present comparison is offered in the hope that it provides clarification and ideas to 
all interested in modeling the coupled thermoelasticity of plates and similar structures. The 
emphasis here is on the mechanical accuracy of the various equations, and in order to compare 
them we consider the most important dynamic feature introduced by the coupling to the 
temperature field: damping of otherwise nondissipative flexural waves. The paper begins 
with an overview of the equations common to all four models in Section ^ Subsequently, in 
Section |21 the models of Chadwick, Lagnese and Lions, Simmonds and Norris are introduced 
in a coherent framework, with some new results for the latter three. Section [3] contains a 
parallel analysis of the dispersion relation for quasi-flexural waves according each model, and 
a detailed discussion and comparison of the damping predicted by each is given in Section 



1 Common equations 

The system of equations derived by Chadwick [T] include the other three models as different 
approximation. In order to see this we first describe what they all have in common. We 
begin with the elastic equilibrium equations with no applied forces 



m 



&33,z + 03/3,0 

0737,7 + cr 3( g j2 



PU3,U, 

pup,tt, /3 = 1,2. 



(la) 
(lb) 



The displacement is u(x, y,z, t) = (ui, u 2 , u 3 ), er is the symmetric stress tensor, p is mass 
density per unit volume, and the repeated Greek suffices indicates summation over 1 and 2. 
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The plate of thickness h occupies — oo < x, y < oo, —h/2 < z < h/2, and is traction free on 
the surfaces: <7i 3 = cr 2 3 = °33 = on z — ±|. 

Integrating eq. (jlaj) with respect to z, using the fact that the shear stresses are constant 
(zero) on the surfaces, and that the normal stresses vanish, implies 

h/2 

dz za 3 p !Z p = phwtt- (2) 

-h/2 

where w is the averaged through-thickness displacement 

h/2 

w(x,y,t) = - J dzu 3 (x,y,z,t) . (3) 

-h/2 

The in-plane equilibrium equations (|lb|) then imply 

h/2 

Mfafr = phw tt + p J dzzufcpa, (4) 

-h/2 

where the first moment of the stresses, or more simply, the moments, are 

h/2 

Mp,{x,y,t)= J dzzapy, (3, 7 = 1,2. (5) 

-h/2 

We now assume the kinematic ansatz normally associated with the name Kirchhoff, 

up = -zw t/3 , (3=1,2, (6) 

so that (JH) becomes 

Mp^py = phwtt - pIAwu, (7) 
where / = h 3 /12, and A is the in-plane Laplacian, A = -J^ + -§^- The development so 
far is independent of the thermoelastic properties of the plate. In fact, equation (J7J) is still 
exact apart from the simplification of the the final term - the rotatory inertia - which is 
often ignored in a first approximation. We will do this later, but for the moment retain this 
inertial term, associated with Rayleigh 

The thermal properties enter through the equations for the in-plane stresses in an isotropic 
thermoelastic plate, which follow from eq. (6.5) of [2] or elsewhere, 

E . aE 

E , x olE . . 

o"22 = n(f£n + £22)-- 0, (8) 

1 — v A 1 — V 

E 

°"i2 — ~ — ; — £12 ■ 
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Here £g 7 = (w/j i7 + w 7i/ g)/2 are the in-plane strains, and #(x, 2, t) is the temperature devia- 
tion from its ambient constant value T. Also, E is Young's modulus, v Poisson's ratio, and 
a is the coefficient of thermal expansion. Combined with eqs. (JSJ) through ((7|), these imply 
the first of our governing equations: 

DA 2 w + Da(l + u)AQ + phw tt - pIAw tt = , (9) 

where D = EI /(l — i^ 2 ) is the bending stiffness, and Q(x,y,t) is the (normalized) first 
moment of temperature: 

h/2 

1 f 

0(x,y,t) = j J dzz9[x,y,z,t) . (10) 

-h/2 

Equation (jH)) is the standard plate equation with the added term which couples the 
plate deflection to the thermal effects. The four theories all use equation in one form 
or another, although an additional equation is required to close the system. The theories 
considered differ in how the temperature moment O is calculated in each. Thus, Simmonds 
[3] and Norris |3] give explicit though different forms for 0. In the theory of Lagnese and 
Lions [2] the function Q(x,y,t) satisfies a 2+1 dimensional partial differential equation. 
Chadwick's model £Q requires solving a 3+1 dimensional equation for 8(x,y, z,t), and as 
such, is the most general of the theories. We begin with this. 



2 The four theories 

2.1. Chadwick's thermoelastic flexural plate equations 

Chadwick [T] derived a coupled system of equations for w(x,y,t) and 9(x,y, z,t). We begin 
with the relation for entropy s in terms of strain e and temperature deviation from ambient, 



1-2^ T 
and the heat flux condition 

Ts t = K(A6 + 6 zz ). (12) 

Here, C and K are the heat capacity and thermal conductivity, respectively. Equations (jllj) 
and ()12|) are well known, e.g. [IB], and are eqs. (6.16) and (6.19) of |2j. Under the thin plate 
assumption, the stress 033 vanishes and therefore 

y ie 11 + e 22 ) + ^9. (13) 



xx T 

Substituting from this into (|TT|) and (|12J). combined with the Kirchhoff assumption (JHJ) yields 

C9 t - K(A9 + 9 ZZ ) - ^^Aw t = 0. (14) 
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This is the second equation of Chadwick, which, together with the deflection equation (J5J), 
forms a closed pair of partial differential equations. The remaining ingredients are boundary 
conditions; however, since we consider the plate of infinite lateral extent, no edge conditions 
are necessary. The PDEs (JHJ) and ()14|) need to be supplemented by conditions for the 
temperature at the surfaces. 

In Chadwick's analysis, and in the later comparisons, the rotatory inertia in (jH)) is ignored. 
In addition we consider the simple no-flux conditions for temperature on the plate surfaces. 
With these clarifications, the Chadwick model may be summarised as 

DA 2 w + Da(l + u)AQ + phw tt = 0, (15a) 
znTE 

C6 t -K(A8 + 6 ZZ )-- Aw t = 0, (15b) 

1 — v 

9 z (x,y,±^,t) = 0. (15c) 

Chadwick's derivation included additional details which gave careful account of whether 
quantities are isotropic or isentropic 1 . We ignore this distinction here but refer readers to 
the original paper for a proper description of the governing equations. 

2.2. The Lagnese and Lions model 

The Lagnese and Lions [2j equations for flexural-thermal dynamics of a thin plate are 

DA 2 w + —a(l + u)AQ + phw tt - pIAw tt = 0, (16a) 

h aTE 

C& t - K AG + -K - — Aw t = 0. (16b) 

These equations follow from, for instance, eq. (6.32) of [2] with zero forcing (p = fs = in 
[2]) and with zero- flux boundary conditions on the top and bottom faces of the plate (Ai = 
in j2]). The first equation ()16a|) is a modification of the standard equation of flexural motion 
for a thin plate , while eq. ()16b|) for the temperature moment, Q(x,y,t), is specialized 
to the case of zero thermal flux on the faces z = ±h/2. Lagnese and Lions [2] included 
the possibility of more general flux conditions. The final term in (|16aj) represents rotatory 
inertia, analogous to Rayleigh's theory, and is a common addition to the classical theory 



The equations as derived [2] contain two errors which need to be addressed before the 
model can be compared with others. First, the factor of 1/2 in the second term of ()16aJ) is 
incorrect, and should be unity. The error is a result of improper application of the variational 
derivation of the equation. The error is evident by comparison of ()16aj) with 0. Appendix 
A provides details on the necessary correction to the variational formulation. The second 
error is the factor (1 — 2v)~ 1 in the final term of (j!6b|) . which should be (1 — ■ The 
error arises from failure to properly apply the entropy constitutive relation to the thin plate 
deformation as done in eqs. (jllj) and (fT3*|) (see eq. 6.19 of |2j). 

1 Chadwick employed a different definition of the thermal expansion coefficient a, i.e. 1/3 of the quantity used here. 
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Equation (j!6b|) may be obtained as the first moment of eq. (fT3j) divided by /. The 
term (h/I)KQ in eq. ()16b|) comes from the evaluation of (K/I)p where p(x,y,t) is the first 
moment of —9 ZZ Applying the assumed zero flux conditions on the plate faces ()15c|) . the 
moment reduces to, with no approximation, 

h/2 

fi = - / dzz9 zz 



-h/2 

= 9(x,y,^t)-9(x,y,-^t). (17) 

In order to relate this to 0, Lagnese and Lions assume that 9(x,y, z,t) is linear in z (even 
though this implies that 9 ZZ = 0), i.e. 

6(x,y,z,t) nzQ(x,y,t), p = hO, (18) 

and hence the term (h/I)KQ in eq. (j!6bj) . 

For the sake of simplicity we choose to ignore the rotatory inertia term pIAw tt in order to 
compare the Lagnese and Lions equations with the other thermoelastic thin plate theories, 
which is not central to the Lagnese and Lions model, and is normally ignored in other first 
order thin plate models. This term could be included but presents algebraic complications 
that make comparisons between models unnecessarily clumsy. With this said, and with the 
corrections noted above, the following system represents the Lagnese and Lions model for 
further discussion: 

DA 2 w + Dct(l + z/)A9 + phw tt = 0, (19a) 
CQ t -KAQ + jKQ-^—Aw t = 0. (19b) 



2.3. Simmonds' model 

Simmonds proposed his model as a simpler alternative to the Lagnese and Lions equations. 
The final result is a single governing equation, which we rederive here using slightly different 
methods. Simmonds' analysis is premised on certain assumed scalings, which for our pur- 
poses are unnecessary as we will retain the leading order terms by physical argument. In 
the process of his derivation, Simmonds calculates the shear stresses in the plate. Since this 
is in itself a useful result, we begin there. Thus, fllb|) and the in-plane stresses of (jSJ), imply 

E 1 1 

03/3,* = ~ 1 _ 2 b(i + u ) u iaP + 2 (- 1 ~ u ) u Pa~i - i 1 + u ) a0 A + P u P±t ■ (20) 

This is comparable to eq. (3.5) of [3], although the latter does not include the final term 
in (|2U|) . which we retain in order to be consistent with the general derivation of Q. This 
term is absent in Simmonds' analysis by virtue of the scalings chosen. Using the Kirchhoff 
assumption of ©, 

osp,z = _ 2 [zAw t/3 + (1 + v)oi9 t p] - pzw^tt , (21) 
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which can be integrated, to yield an explicit equation for the shear stresses in terms of w 
and 9, 

1 / 9 / E . Ecu /",/„/ /x /x 

^ = ^ 2 -—)(^-—^Aw t p-pw tm ) + ^—^ J dz'9^x,y,z',t). (22) 

-h/2 

This clearly satisfies the conditions that the shear stresses vanish on the surfaces as long as 
9(x,y,z,t) is an odd function of z, which is assumed. Equation is analogous to (4.11) 
of [3] although the latter contains an explicit form for the integral of 9. We will return to 
this point in detail later, and leave (}2*2*j) in its general form at this stage. Substituting from 
(EH) gives 

DA 2 w - pIAwu- J dz J dz ' & d ( x ,y, z ', t )+ P hw tt = 0, (23) 

-h/2 -h/2 

and integrating by parts reproduces eq. (jHJ) exactly. The Simmonds model therefore in- 
cludes the basic equation for the plate deflection. It remains to discuss his equation for the 
temperature. 

Without going into the details of the scaling employed, it suffices to note that Simmonds 
solves eqs. (|15b|) and ()15c|) with A9 — > 0, i.e. 

ZZ - %6 t + Aw t = 0, (24a) 

K Kyi — v) 

9 z (x,y,±^,t) = 0. (24b) 

These are equations (4.3) and (4.4) of Simmonds' solution may be expressed in a slightly 
different form than the original (eq. (4.9) of 0). Consider the ansatz 

00 1z 

9(x, y, z,t) = Yl mW ( T ) v {n) (x, y, t) + 9 H , (25) 

n=l n 

where the series represents the particular solution, and 9h is a solution of the homogeneous 
equations which is required in order to give the correct initial condition for the temperature. 
The latter may be expressed in terms of eigenfunctions, as in 0. The functions u {n \z) are 
defined such that: 

u (")"(s)=« (n - 1) (s), u (n) '(±l)=0, n>l; u^(s) = s. (26) 

These may be obtained recursively as described in Appendix B. Note that 9 of (|23j) auto- 
matically satisfies the zero flux conditions on the surfaces. 

Substituting from (j2*3j) into ()24aj) and equating the coefficients of to zero implies 
recursion relations for v^ n \ 

v^ +1 \x,y,t)= f ^4 n \x,y,t),n>l, „«(*, y, t) = . (27) 
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Thus, formally at least, 

f^^-^^i)^" 4 ' + °" (28) 

where 

- - % < 29 > 

is a characteristic time for thermal diffusion across the plate thickness 2 . The first moment 
of can be evaluated using eq. fl2SJ), with the result 



Aw; + Q H , (30) 



O T"l p CO Q 



oTE 



2d 17 d 2 62 9 2 
5 T5 9t ~ 105^9? + 945^ 5 dr 



with obvious meaning for 0#, which we assume is zero for simplicity. The coefficients 
follow from the general formula (|B.6|) in Appendix B. The convergence properties 
of the series expansion are also discussed in Appendix B. 

Simmonds argues that the times of interest are much longer than t$- Therefore only 
the first term in the expansion of the particular solution is significant, yielding an explicit 
expression for the first moment, 

h 2 aTE 

B = 77-— rAu> t , bimmonds. (31) 

lUK [I — v) 

This implies a single governing equation for the plate deflection (eq. (4.13) of 0) 

DA 2 w + e u ^T S DA 2 w t + phw tt = 0, (32) 
5 

where the nondimensional factor e u is defined as 

eu = —p^(z ■ 33 

C 1 — v 

We will use this quantity repeatedly in the remainder of the paper. 

It is interesting to note that of all the theories Simmonds' does not involve the heat 
capacity C in the final equation. This is a consequence of the fact that the approximation 
(|33j) is equivalent to ignoring the time derivative of 9 in ()24a|) and solving instead 

zn/TE 

K9 ZZ = -- Aw t , (34) 

1 — 1/ 



subject to the zero- flux conditions (|15c|) . The solution of this is obviously independent of 
the heat capacity, and is 

z ,h 2 z 2 , aTE . , , 

= -(---)— —K Wt . (35) 



2 Note that rg = tq, where to, defined in eq. 1391 . is the fundamental time scale for all the models. 
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Substituting this approximation for the temperature distribution into (JTUJ) yields (pTTj) di- 
rectly. 

2.4. Norris' thermoelastic flexural plate equation 

Norris' model is a single viscoelastic-type equation for plate deflection w(x,y,t), 

(1 + e v ) DA 2 w -e u Dg* A 2 w + phw tt = 0, (36) 

where * denotes convolution in time and e„ is defined in (133)1 . The relaxation function g(t) 
is 

9(t) = j:^e-^, (37) 

with 

= 7r 4 (2n + l) 4 ' ' n_ (2n + l) 2 



n=0 Tn 



d n = 4/0 96 , r n = T ° , n = 0, 1, 2, . . . , (38) 



and 

-o = % (39) 

The single equation (jSfij) may be obtained in a variety of ways. For instance, it is clear 
from Appendix B that it is equivalent to retaining the complete series in eq. (130)1 rather 
than just the first term as Simmonds does. A more direct derivation is outlined in Norris 
and Photiadis [H] in which the time-harmonic version of (J3fi)l is derived. Briefly, this follows 
by solving ()24a)l in the frequency domain and converting back to the time domain (see eqs. 
(jfgj) and dHH) yielding 

olTE 

= eft _ u ) ( Aw ~9* Au>) . (40) 

When combined with (jl5aj) . this reduces the plate model to eq. (136)) alone. The justification 
for ignoring the in-plane spatial derivatives in (j!5b)) is motivated by clear physical arguments 
in [H] and is also evident from the results below where we will see that predictions of the 
Chadwick and Norris models are essentially identical. See [3] for a more thorough and 
complete derivation of the simplified thermoelastic plate theory, along with appropriate 
plate boundary conditions and example applications. 

2.5. Summary of the four models 

All models share eq. Q, which we simplify for present purposes by ignoring the rotatory 
inertia, so that 

DA 2 w + Da(l + u)AQ + phw tt = 0. (41) 
The theories differ in how they estimate of ()10p . Thus, 

C6 t - K{A9 + 6 ZZ ) - Aw t = 0, with 6 z (x, y, ±-, t) = 0, Chadwick, (42a) 
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7 rji 771 

C& t - KAQ + -KB - Aw t = 0, Lagnese & Lions, (42b) 

I 1 — v 

2 

jfjToAwt, Simmonds, 
9 = . - . x { (42c) 



a 1 + v) 



Norris. 



In Chadwick's model B is determined indirectly by first solving the equation for 8(x, y, z, t), 
while 0(x, y, t) is found directly as the solution of a PDE in Lagnese and Lions' model. The 
Simmonds and Norris models provide explicit expressions for 0, although in the latter case 
a convolution is required. We now look in detail at solutions to these equations. 



3 Flexural waves and attenuation 



3.1. Traveling wave solutions 

In order to compare the four models we disregard any boundary effects and consider the 
infinite plate. The appropriate fundamental solutions are traveling time harmonic straight 
crested waves. We examine the form of these waves for the different theories with a view 
towards distinguishing the models and determining a best approximation to the underlying 
physics. The assumed form of the solution is 



{w(x,y,t), 6(x,y,z,t)} 
{w(x,y,t), Q(x,y,t)} 
w(x,y,t) 



{w , 0ou(*)}e <(feB - wt) , Chadwick, 
{w , e }e i(fcc ~ a; * ) ) Lagnese & Lions, 
Wo e i(kx-ujt) ^ Simmonds, Norris. 



(43a) 
(43b) 
(43c) 



The frequency u is assumed to be real and positive. Define the real- valued flexural wavenum- 
ber k > and the complex-valued thermal wavenumber 7 by 



4 2P^ 2 

k =u —, 7 



IUO- 



c_ 

K' 



(44) 



The dispersion relation between k and u has a common form for in all four plate theories 



k 4 



1 + 



K 



0. 



(45) 



where F is a different quantity for each, 

f F c (( 7 2 - k 2 f/ 2 h), Chadwick, 

Fl((i 2 — k 2 y/ 2 h), Lagnese and Lions, 

Simmonds, 

, Fc(jh), Norris, 



10 ; 



(46) 
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and the functions F c and F L are 





x 2 - 12' 



1 



(47) 



These results follow from the identity jUj 



oo 




f(jh) = 1 - g(u), where = 

n=0 



(48) 



and the g{uf) is the Fourier transform of g(t). 

All of the models considered are deviations from the non-dissipative elastic Kirchhoff 
model of thin plate dynamics. Any changes in wave speed, energy propagation speed, real 
wave number, modal frequency, or a quantity associated with energy preserving systems, is 
expected to be a perturbation of the value without thermal effects. It makes sense to consider 
and compare the effects of the different models on the attenuation as this is an order unity 
effect. While there are several measures of dissipation, we focus here on the quality factor, 
Q, a non-dimensional parameter, defined in the present context by the imaginary part of the 
wavenumber as 



The quality factor can be related to, for instance the attenuation per unit length of a prop- 
agating wave as uo / (2vQ) , where v is the speed of energy propagation. 

3.2. Perturbation expansion 

The dispersion relation can be expressed as a cubic in k 2 for the Lagnese and Lions model. 
However, rather than solve this we take advantage of the fact that for materials of interest 
e„ < 1, see Table I. This permits solution of (|45|) by regular perturbation. The leading order 
roots correspond to the classical theory, i.e. k 4 — k 4 = and we consider the propagating root 
with unperturbed solution k = k. The leading order solution is therefore the non-dissipative 
flexural wave with real phase speed u/k and real energy propagation speed 2u/k. The first 
order term in the perturbed solution is given by 




(49) 




(50) 




K 2 y/ 2 h), Chadwick, 



K 2 ) l l 2 h), Lagnese and Lions, 



(51) 



10' 



Simmonds 



Fc(lh), 



Norris. 
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Material 


e v x 10 3 


I 

(nm) 


Aluminum 


8.38 


32.3 


Beryllium 

.. 


5.91 


11.2 


Copper 


5.71 


91.8 


Diamond 


0.13 


58.1 


GaAs 


0.96 


230 


GaPh 


0.69 


38.1 


Germanium 


1.08 


25.7 


Gold 


0.42 


143 


Lead 


19.11 


54.1 


Nickel 


5.15 


13.6 


Platinum 


2.92 


24.9 


Silicon 


0.40 


28.1 


SiC 


0.87 


5.81 


Silver 


7.50 


177 



Table I. The thermoelastic coupling parameter of eq. (|3*3*|) and the phonon mean free path / 

of eq. USSJx, both at at 300 K. 

The replacement k — > k in (j51|) is important as it leads to further significant simplification, 
as we see next. Specifically, eq. ()50|) may be replaced by 



?/(7h)+0( £ ;,£„«), 



*=1- 

re 4 

where 5 <C 1 is defined below, and the function / is 



(52) 



/(*) 



1 + 24 f£ 



tan §), 



x 2 -12' 



First note that 



x~ 
10' 



— = -iaS, 

T 



Chadwick, Norris, 
Lagnese and Lions 
Simmonds. 

where S = ~, 



(53) 



and the length / and nondimensional parameter a are defined as 

2 r /— r (1-u) 



cC 



2 1 



(54) 



(55) 



Here, c is the averaged elastic wave speed, c — (cl + 2ct)/3, where = £'/[2(l + u)p], 
c 2 L = E(l — u) /[(l + u)(l — 2v)p\. The length / is defined in this way so that it is identifiable 
as a thermal phonon mean free path, since c is a typical phonon speed jlTj. The parameter 
a is of order unity. Values of / for a variety of materials are given in Table 1 . The phonon 
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mean free path must be far less than a typical structural length scale if thermoelasticity 
theory is to remain valid, and in particular we must have I <C h. Consequently, 

S < 1, (56) 

and 

( 7 2 - K 2 fl 2 h = yh [1 + 0(5)] . (57) 

This provides the justification for the simplification (J52j) and (|53|). Alternative arguments, 
not based on the phonon mean free path, but which give the same result are given in 
Appendix C. 



3.3. The attenuation 

The leading order approximation to the damping of the thermoelastic plate wave is, based 
upon pnjl . (JHUJ) and the approximation for F , 

Q- x = -e„lmf(>yh). (58) 

The quantity /(7/1) is simple for the Lagnese & Lions model and Simmonds' model, for 
which (|58jl is exactly 



Lagnese and Lions, 



Q- 1 = e„ x ^ (59) 
w r o Tn > S immonds , 



where 



7T 2 



7 i=12 7 b. (60) 

Turning to the damping as predicted by Chadwick's model, which is the same as for 
Norris' by eq. (JS3J), eqs. (|58|) and (f4"%|) imply 

00 ^ 

Q~ x = e v - — - — 5- 9) Chadwick, Norris. (61) 

n=0 1 + U T n 



4 Comparing the four models 

4.1. Comparison of attenuation 

In order to compare the formulae for the damping, define 

p _ d n UJT n LUT L IT 2 UT 

h^^r L ~TT^vT ~ UTo To' s °-iT^ 



Thus, Q~ x is equal to e u S for the Chadwick and Norris models, e u Si for Lagnese & Lions, 
and e v Ss for Simmonds' model, respectively The single term expression So will be explained 
presently 
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Of the four predictions, Simmonds' is the only one which gives damping that grows 
without bound as u — > oo. The Simmonds equation is not of the form of a standard 
relaxation, but more like a pure viscous damping effect. We note however, that at low 

oo 

frequency, S = Ss+O((ujt ) 3 ), which follows from the identity d n T n = r 7r 2 /10. The 

n=0 

model predicts the right linear dependence as ujt — ► 0. In summary, Simmonds's model is 
correct at low frequency only. 




4 



Figure 1: The functions S, Sl, and So of eq. H(>'2l . The right hand frame focuses on the peaks near wro = 1. 



The three expressions S, Sl, and So are compared in Figure [T] as functions of ur . It is 
evident from Figure Q that while the Lagnese and Lions function Sl has the same overall 
behavior as the more precise prediction S, the latter is better approximated by Sq. The 
maximum absolute error in approximating the sum S by Sq is less than 6 x 10~ 3 . By 
comparison, the absolute error obtained by approximating S with Sl is on the order of 10 
times larger. More serious than this, perhaps, is the fact that the frequency at which the 
maximum of Sl occurs is noticeably in error, at utq = 12/ n 2 = 1.216, see Fig. ^ 

It is interesting to note that the maximum of S occurs close to but not exactly at ur = 1. 
The peak occurs when 

S IT 2 

= _ = 0.4935. (63) 

lotq 20 V ; 

This must be solved numerically for u;r , which can be done using with (J62j) or the following 
alternative formula for S, which follows from (|53jl . and from [33], 



x 6 



x — 



f sin x + sinh x 
Vcosx + coshx 



7T 



(64) 



We find that the maximum is at lutq = 1.0014143 . . ., and the corresponding maximum value 
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is "S'max = 0.494178 . . .. The approximant Sq has its maximum precisely at ujt = 1, and its 
maximum value (1/2) is slightly in excess of the correct value S m3X by about 1%. 

4.2. Discussion 

We can gain further understanding of the above results by returning to the governing plate 
equations. The approximation (|57|) may be interpreted in terms of the in-plane spatial 
derivatives in both (J19bj) and Ijlfibj) . Any variations with wavenumber on the order of k or 
longer can be safely ignored because the thermal diffusion makes the variation in z dominant. 
Therefore, for flexural waves and most disturbances associated with thin plate dynamics, we 
can ignore the A0 in (j!9b|) and AO in (|15b|) . With the latter approximation, eqs. (j!5b|) 
and (|15c|) can be solved to yield eq. (|40|). as discussed before. In this way the Norris 
model and the Chadwick model give essentially identical predictions. It also indicates that 
the Chadwick plate model, while correct, is unnecessarily cumbersome: the in-plane spatial 
derivatives of temperature can be safely ignored jH]. By the same argument, we may ignore 
the term KAQ in (j!9b|) . so that the correspondingly simplified form of the Lagnese and 
Lions model is 

CQ t + -KQ = ^ Aw t . (65) 
n z 1 — v 

This can be solved as 

aTE 

= Tvl V ( Auj " 9l * Aw) , (66) 

Ly{L — V) 



where 



g L (t) = -e- t ^. (67) 



Substituting from ()66|) into ()19a|) it follows that the Lagnese and Lions thermoelastic plate 
model reduces to 

(1 + e„) DA 2 w -e u Dg L * A 2 w + phw tt = 0. (68) 

Referring to eqs. (|37|) and (|62j). it is clear that the one term damping function Sl can 
be associated with the single relaxation function ^(t) of (|67|) . At the same time, the 
approximation S — > So is equivalent to replacing g(t) of (|37|) by a simplified relaxation 
function with a single relaxation time, g(t) —>■ go(t), where 

9o(t) = -e*to. (69) 

Therefore, rather than comparing the frequency domain functions, S, Sq and Sl, one can 
equally well look at the time dependent relaxation functions g(t) of ([37)1 . gL(t) of ()67|). and 
9o{t) of These are compared in Figure El 

The relaxation functions g, go and gL share the common property that their integral is 
unity, 

oo oo 

J g {t)dt = J g L (t)dt = 1, (70) 





J g(t)dt 
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g (t) 



Figure 2: The function rog(t) compared with Togo(t) and TogL{t) 



but the three functions have distinct first moments, 

oo 2 00 

J tg(t)dt= 7 ^T , J tg (t)dt = r , 



oo 2 

y tg L {t)dt = ^r . 



(71) 



The functions and giif) have the same values at t = 0, g(0) = «7l(0) = 12/(7r 2 r ) while 
«7o(0) = 1/ T 0) as evident from Figure 121 The identities (|7U|) and (|7T|) are associated with 
the first two terms in the low frequency expansions of the Fourier transforms. The exact 
transforms are g, defined in (JUL and go and cjl, given by 



1 



9l{uj) 



1 — ibJTQ ' 

and the low frequency expansions are 
g(u) = 1 + iuj r — + ... , g (uj) = 1 + iujtq + ... 



1 — IUJTl 



71 

°12 



(72) 



(73) 



Hence, one reason that ^o(^) an d >S'o( a; ) are much better approximations than giS$) and Sl(u) 
to the precise quantities g{t) and S(a;), respectively, is that the number n 2 /10 = 0.987 is 
closer to unity than 7r 2 /12 = 0.822. 

These results imply that the Lagnese and Lions model could be improved by replacing 
the term (h/I)KQ = (12/ h 2 )KG in ijTTSKjl with (10/ h 2 )KQ. The resulting simplified single 
plate equation, analogous to Norris's equation (J36j) is then 



DA 2 w -e u Dg * A 2 w 



phwu = 0. 



(74) 
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Note the relaxation function is now g , rather than g L as in eq. (|bS|). The use of produces 
the inaccuracies shown in Figures ^ and HI while the same Figures illustrate that go is a f ar 
better approximation to g(t). 

The derivation of eq. ()74|) as an improvement on (|UH|) indicates that the source of the 
discrepancy in the Lagnese and Lions model can be partly ascribed to the term (h/I)KQ 
in eq. (|19b)l . As noted before, the term in question originates in the evaluation of (K/I)n 
where the first moment of — 9 ZZ is defined in (|17|) and evaluated, according to the linear 
approximation for 9 assumed by Lagnese and Lions, in (|18|). Equation (J18)) is incorrect. The 
correct form follows from Norris and Photiadis [H], and is most readily expressed in the 
frequency domain, as 

~ 1 sin *yz ~ 

9(x,y,z,u) = - ^—r{z K )Q(x,y,u), (75) 

l-g{uj) 7cos7f 

where 7 and g{uj) are defined in eqs. (|4*4"|) and (jIBJ), respectively. Thus, eqs. |T7]) and (f7H|) 
imply 

and expanding, 

5 _ (jKf_ 4 4 h ~ 

[ 6 1008 U ;J V ; 

This suggests that a better approximation than that of Lagnese and Lions is to take the 
leading order term, 

/i = ^e, (78) 

o 

which in turn implies the modified governing equation (J74J) . 

What are the relative benefits of the four different models? The results presented indicate 
that for the infinite plate, one might just as well use the Norris model, since it is almost the 
simplest (only Simmond's is simpler) yet it provides the same accuracy as that of Chadwick. 
The equivalence rests on the fact that for any practical plate of thickness h, one must have 
h 3> I, the mean free path of the thermal phonons. It might be argued that the other 
models are more appropriate for finite systems, i.e. plates with edges subject to a variety 
of boundary conditions. We have not discussed the edge conditions for each model, as this 
is a large topic in itself [21 II]. Suffice to say that for each model one can define appropriate 
conditions that supplement the governing equations and define a closed problem for a given 
practical situation. 

In order to examine the utility of the models for a finite plate, Appendix D discusses 
the solution for the case of a beam with fixed ends and zero heat flux at the ends within 
the context of the Lagnese and Lions model. The results there indicate that the additional 
equation (jl9b|) for B only effects the solution in a boundary layer at the edges. The Norris 
and Simmonds models do not provide a boundary layer solution, but have solutions analogous 
to the classical Kirchhoff theory but with modified (complex valued) flexural wavenumber. 
Once could argue that the Lagnese and Lions model is useful in that it provides a boundary 
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layer, even though the layer, being exponentially narrow, has insignificant effect on the 
value of both the modal frequencies and damping of the oscillating modes (see for further 
details). The presence of a boundary layer at the edges is expected, and while the Lagnese 
and Lions model predicts one, only the Chadwick model is capable of providing the proper 
description of the boundary layer. Thus, the width of the boundary layer of the Lagnese 
and Lions model is 0(S 1 ^ 2 /k) where 5 is defined in eq. (|54|) and k is the classical flexural 
wavenumber. This boundary layer width could be less than, on the order of, or larger than the 
plate thickness, depending on the material parameters. One can envisage circumstances in 
which the boundary layer width is less than h, and consequently will have z— dependence. 
This is not contained within the Lagnese and Lions model and can only be captured by 
the Chadwick model. This fundamental difficulty arises with all beam and shell theories, 
classical or higher-order: they do not accurately describe edge boundary layers that vary 
over a distance 0(h). However, regardless of the precise nature of the temperature field 
near the boundary, it will have essentially no effect upon the modal properties of the plate 
(frequency, damping, mode shape, etc.). 

5 Conclusion 

We have compared four distinct models for flexural motion of thermoelastic thin plates, 
with specific attention to the leading order effects on wave damping. Of the four, the 
Chadwick theory is expected to be the most reliable as it retains the 3D nature of the 
thermal field. While the Lagnese and Lions has the same overall behavior as Chadwick's 
in terms of frequency, there is a distinct quantitative difference. This is attributable to the 
approximation of the first moment of the temperature, which is inaccurately accounted for 
by the Lagnese and Lions model. The relatively simple single equation of Simmonds predicts 
the correct low frequency behavior but is incorrect for medium and high frequencies. The 
reason is that the Simmonds equations approximates the thermoelastic effect as a simple 
viscous damping term which does not capture the full frequency content. Comparisons of 
the Norris model with that of Chadwick indicates that the former adequately predicts the 
full frequency behavior. The single equation of Norris provides an accurate approximation 
because the the dominant variation in the temperature field is in the through thickness 
direction, which permits it to be replaced by an equivalent viscoelastic term. The Norris 
model retains the simplicity of a single uncoupled equation for plate deflection, but provides 
a more accurate representation of the thermoelastic "hidden variables" than Simmonds'. 
The Lagnese and Lions model misses the crucial point that the dependence of 9 on the in- 
plane coordinates x and y is secondary to its dependence upon z. This permits considerable 
simplification with the result that a separate equation for temperature is not necessary. 
Instead, the temperature is incorporated into the classical plate equation as an effective 
viscoelastic relaxation [3] . Of the four models, the simpler single equation of Norris captures 
the thermal mechanics of the more general Chadwick theory, and is to be preferred to the 
other models considered. 
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Appendix A: A correction to the Lagnese and Lions equations 

In this Appendix the correction to eqs. (j!6a)l is obtained by the same method used in the 
original derivation of Lagnese and Lions [2] . The cause of the Lagnese and Lions error is that 
the energy density U (different notation) is defined in [2] by U — \cr ■ e. The correct energy 
density for the conjugate independent variables strain and temperature is the Helmholtz free 
energy F, see e.g. JEJ, where F — |(<x • e — s9). Using equations (fT3|) and (JBJ), the total 
strain energy becomes 

F=\j [d[(Aw) 2 - 2(1 - v)(w n w 22 - w 2 12 ) + 2(1 + v)a9Aw] - -6 2 }dxdy . (A.l) 
n 

The analogous equation (6.7) in j2j omits the final term and the coefficient of the (l + v)aOAw 
term is unity rather than 2, and hence the error in (jlfiaj) . 



Appendix B: 

The functions defined by (|2fi|) may be found recursively as follows, 



) = Va"— ■ , (B.l) 



2n-2j+l 



2n -2j + l 
with 

qU— 1 n— 1 

a? = — ^rA — — , n > j > 0; < = - ^ a?, n > 0; a° Q = 1. (B.2) 



(2n-2j)(2n-2j-i; 



j=0 



The values it^(l) in eq. (|3()|) can be found using these recursion relations. We now 
demonstrate another, more constructive, way to obtain these coefficients. The series expan- 
sion (|3U|) may be obtained from the Norris solution for written in the frequency domain 

as 

aTE 

6(x, y, u) = c ^_ v) [1 - 9{u)\ A«5, (B.3) 
where g(u) is given in (j48jl . The latter may be expanded in a Taylor series as 

oo oo 

Sf(u;) =J2 S m where S m = £ rf n r n m . (B.4) 

m=0 n=0 

The infinite series 5" m may be evaluated using (|38|). and known results j3H], as 

C ^™ /IO„2m/r)2m+4 ,\ |-B2m+4| /T3 _\ 

^ = T ° 487F (2 _1) (2^T^!' (R5) 
where £> m are the Bernoulli numbers jlH]- Eliminating To in favor of Ts of (|29|) gives 

^ m -3r 5 2 (2 -1) (2m + 4)!' [ ' 
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and hence 

2 IT 62 

fif(w) = 1 + g iwr s + — (iwr s ) 2 + — (iwr s ) 3 + . . . . (B.7) 

Equation (|3T)J) then follows from ()B.3J) . ()B.7jl and the connection — zcj — > <9/<9t. Finally, 
convergence properties of the series (|B.3|) i follow from the behavior of the Bernoulli numbers 
for large index. Omitting the details, it may be shown that the series is unconditionally 
convergent for ujtq < e. 



Appendix C: Alternative justification of the Fq approximation 

The justification for the approximation of F$ (eqs. (|51|) and ([53))) is discussed further. If the 
reader is uncomfortable with the argument based on the requirement that the phonon mean 
free path must be much shorter than the plate thickness, an alternative reason is given here. 
We focus directly on the damping defined in (|49|) as given by eq. (!5(Jj) . Thus, 

Q" 1 = -€ u 7i 2 ur ReF , (C.l) 

where F is given by (|51J1. 

The Lagnese and Lions model is considered first. It follows from (|C.1|) that 



Q 1 = e u ^— r-, Lagnese and Lions, (C.2) 

1 + UJ Z Tf 



where 



7r 2 r 



12 + {nhf (C - 3) 

The thin plate model, whether thermoelastic or classical, is premised on the long wavelength 
assumption that 

nh < 1. (C.4) 

Therefore, 

n = ^r + O[(Kh) 2 ], (C.5) 

and it is entirely consistent with the thin plate theory (Kirchhoff approximation) to take 
r* = r 7r 2 /i2, i.e., r* = T L . 

For the Chadwick model, eqs. (|4~9|). (|5*njl and (|4"5|l imply 



n=0 ± + U T n* 

with 

T, 



Q _ x = ^ ^* Chadwick, (C.6) 
„ „ 1 + cc^r£. 



1 _|_ 2St 

TO 



Invoking the necessary condition (jC.4|) . we have 

T n , = T n + 0(K 2 h 2 ) , (C.l 
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and consequently, to the same order of accuracy, i.e. ignoring 0(n 2 h 2 ), 

dn.UJTn 



Chadwick. 



V = 2^ l T 1 ■ 

n=0 1 ^ ^ 'r. 

Appendix D: Solutions for a finite plate or beam 



(C.9) 



We consider a plate of finite extent in the x— direction, — L < x < L, — oo < y < oo, with 
motion independent of y. The ends are assumed, for simplicity, as fixed, with zero heat flux. 
The problem is completely analogous to that of a beam, which is simply the limit of a plate 
of vanishing width in the y— direction jJUJ. The beam solution is obtained by setting v — > 
but for consistency we will stay with the plate problem here. 

For simplicity, we only seek solutions that are symmetric in x, i.e. of the form 



w (x 



t) — Aj cos(kiX — ut), 6(x, t) — piAi cos(/cj£ — out), 



(D.1) 



i=l 



i=l 



where ±fcj, i = 1, 2, 3 are the three independent roots of ()45|) and pi follow from 

kf - k a £vl 2 k 2 



a« as 



The end conditions 
imply 



Pi k 2 a{l + v) a (l + u)(kf- 1 2 +§)- 
w(±L, t) = w x (±L, t) = 9 X (±L, t) = 0, 

cos k\L cos k 2 L cos k%L 
ki sin h\L k 2 sin k 2 L /c 3 sin k^L 



(D.2) 
(D.3) 





-A 1 - 






A 2 


= 0, 




. A 3 _ 





(D.4) 

pikisinkiL p 2 k 2 sin k 2 L p^smk^L 
and in order that we have a non-trivial solution, (u, k) must satisfy the dispersion relation 

cos k\L cos k 2 L cos k^L 



k± sin k\L k 2 sin k 2 L k% sin k^L 
PikisinkiL p 2 k 2 sink 2 L p 3 k 3 sink 3 L 



0. 



(D.5) 



Let k\ and k 2 correspond to the roots which are close to (within order e u ) the elastic roots 
±n and ±in (i.e. the roots of k A — k a = 0). The wavenumber /c 3 is, to leading order, given 
by 3 

(D.6) 



k 2 



7 ~h?- 



3 This follows, for instance, by equating the two expressions for pz in ID. 21 . Note that it does not assume any scaling between 
7 and 1/h, since yh =0(nh/ '8 1 / 2 ), where both nh and 5 are small quantities. 
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It follows that 



where 5 <C 1 is defined in (|54|) . Thus, 



Pi 




P2 






P3 


1 


P3 




mi 2 - m 



\pl\, \P2\ < |P3| 



which enables approximation of eq. (jD.5|) as 

cos/ciL cos k 2 L 

kisinkiL fc 2 sinA; 2 L 



0. 



(D.7) 



(D.8) 



(D.9) 



To be more precise, k\, k-i are roots of k 4 — k' 4 = where the complex valued wavenumber 
k' follows from (|52j). 

K' 4 = K 4 [l-e„f(jh)}, (D.10) 

and / is defined in (|53|). Thus, to leading order in the thermoelasticity the dispersion relation 
is, from (|TX9|l . 

tan k'L + tanh k'L = . (D.ll) 

It is important to note that this dispersion relation has the same form as that of the elastic 
plate or beam without thermoelasticity. The latter is obtained by the reduction k' — > k. 
Thermoelasticity has the effect of modifying the wavenumber in the same way that it is 
altered in the infinite system, and the end conditions do not introduce any new effects. In 
fact, it can be easily checked that the dispersion relation (jD.ll)) is also obtained if one uses 
the Norris equation, which does not include temperature explicitly (and therefore does not 
require an end condition associated with temperature). 

The above result can be explained in terms of boundary layers at the two ends. In order 
to see this from the full solution for the Lagnese & Lions model, we note that p\ and p 2 
may be approximated using (|56|). and combined with ()D.4|) the leading order solution can 
be shown to be 



w 



e 



A ( sinh kL cos kx 



€ u 2k 2 j 2 A 
a(l + ^)(7 2 -# 



3„,2 



sin kL cosh kx 



12' 

7?. 



sin kL sinh kL 



cos k 3 x 
sin k-xL 



Ik 



(D.12) 

cos k 3 x " 



— sinh kL cos kx + sin kL cosh kx + — sin kL sinh kL 

ks sm /c 3 L 



With no loss in generality, choose k$ as the root with positive imaginary part. Then \e iksL \ < 
|e~* fc3i |, and consequently, for positions not near the center of the plate, we have 



e 



€ u 2k 2 j 2 A c 



sinh kL cos kx + sin kL cosh kx 



i ^- sin kL sinh kL e iks (i 1 x 1 } 



(D.13) 
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with a similar equation for w. While (jD.13|) is not valid at the center (x = 0) it provides a 
faithful leading order approximation near the ends. The decay is far more rapid than that 
of the evanescent flexural wave, since lm(ks/K) ^> 1 (see eqs. ()54|) and (|56|)). This indicates 
that the contribution of the complex wavenumber k% is an exponentially decaying boundary 
layer. The fact that \k3\J1 is on the order of or larger than unity raises doubts about the 
utility of the model, particularly near boundaries. 

Acknowledgment: Comments from Prof. James G. Simmonds are appreciated. 
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